#include "interp.hpp"
#include <iostream>
#include <fstream>
#include <vector>
#include <cmath>
#include <functional>
#include <algorithm>

const double PI = acos(-1.0); // Define PI as the arccosine of -1

/**
 * @brief Function f defines the mathematical function 1 / (1 + 25 * x^2).
 * 
 * @param x The input value to evaluate the function.
 * @return The result of evaluating f(x).
 */
double f(double x) {
    return 1.0 / (1.0 + 25 * x * x);
}

/**
 * @brief Computes the Chebyshev zeros, which are the roots of the Chebyshev polynomials.
 * 
 * @param n The degree of the Chebyshev polynomial.
 * @return A vector of Chebyshev zeros within the interval [-1, 1].
 */
std::vector<double> chebyshevZeros(int n) {
    std::vector<double> zeros(n);
    for (int k = 0; k < n; ++k) {
        zeros[k] = cos((2 * k + 1) * PI / (2 * n));
    }
    return zeros;
}

int main() {
    // Define values for n, representing the number of Chebyshev nodes for interpolation
    std::vector<int> ns = {5, 10, 15, 20};
    
    // Perform interpolation for each specified n
    for (int n : ns) {
        // Get Chebyshev zeros (nodes) for the current n value
        auto zeros = chebyshevZeros(n);
        std::vector<double> xValues, yValues;

        // Calculate function values at Chebyshev nodes
        for (double zero : zeros) {
            xValues.push_back(zero);
            yValues.push_back(f(zero));
        }

        // Create a Newton interpolation object with x and y values
        NewtonInterp interp(xValues, yValues);

        // Define output directory and filename
        const std::string output_dir = "output";
        std::string filename = "C_" + std::to_string(n) + ".txt";
        std::ofstream outfile(output_dir + "/" + filename);
        if (!outfile.is_open()) {
            std::cerr << "Failed to open " << filename << " for writing." << std::endl;
            return 1;
        }

        // Evaluate interpolation at points from -1 to 1 with a step of 0.002
        for (double x = -1; x <= 1; x += 0.002) {
            outfile << x << " " << interp.evaluate(x) << std::endl;
        }
    }

    return 0;
}
